Interplay between network structure and self-organized criticality 
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We investigate, by numerical simulations, how the avalanche dynamics of the Bak-Tang- 
Wiesenfeld (BTW) sandpile model can induce emergence of scale-free (SF) networks and how this 
emerging structure affects dynamics of the system. We also discuss how the observed phenomenon 
can be used to explain evolution of scientific collaboration. 
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Since the discovery of self-organized criticality (SOC) 
by Bak, Tang and Wisenfeld (BTW) , the phenomenon 
has received enormous attention among the researchers. 
During these almost twenty years, dozens of variants of 
the original sandpile model 0, were studied 0. [J 0. Q 
and a number of examples of SOC in real world were dis- 
covered. One of the most remarkable features that char- 
acterizes self-organized criticality is the power law distri- 
bution of the characteristic events. The feature combined 
with the abundance of real- world networks with scale- free 
(SF) degree distribution may give rise to sus- 

picion that there exists mutual relation between the two 
issues. Although the idea has been already mentioned in 
several papers (see for example 0,0), according to our 
knowledge, so far no-one has established a link between 
the two phenomena. 

At the beginning let us recap the rules of sand- 
pile model which is a simple intuitive example of self- 
organized criticality. It is a cellular automaton whose 
configuration is determined by the integer variable C, (the 
height of the " sand column" ) at every node i of the net- 
work. Depending on network structure minor differences 
in definition can occur. Here we follow the definition of 
BTW model for random networks with a given degree 
distribution p{k) 0|. The dynamics is defined by the 
following simple rules: A grain of sand is added at a 
randomly selected node i: Ci — ► c,* + 1. A sand column 
with a height Cj > ki , where fcj is equal to degree of the 
node i, becomes unstable and collapses by distributing 
one grain of sand to each of it's ki neighbors. This may 
cause some of them to become unstable and collapse at 
the next time step. This in turn can lead to an avalanche 
of next instabilities. During the evolution a small frac- 
tion / of grains is lost, which prevents system to become 
overloaded. When avalanche dies another grain of sand 
is added. 

In the sandpile model distributions of avalanche sizes 
(measured as total number of topplings in the avalanche), 
avalanche areas (the number of distinct nodes participat- 
ing in a given avalanche), avalanche durations as well as 
many other statistics follow power law distributions. 

Studies of sandpile dynamics carried out so far show 
that the characteristic exponents of measured distribu- 
tions depend on the network topology. For example the 
avalanche size exponent is r = 1 for 2D square lattice 



1] and t = 1.33 for 3D cubic lattice 0]. In Erdos-Renyi 
(ER) random networks r = 1.5 13] . Recently, Goh et 
al. have studied sandpile dynamics on scale free net- 
works p(k) ~ fc -7 |14lll5| and they have shown that the 
avalanche area exponent r is independent on the average 
network connectivity (k) and changes with the exponent 
7 of the degree distribution. They have obtained 
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in the range 2 < 7 < 3 and r = 1.5 for 7 > 3. 
The question we ask in the paper is the following: how 
the avalanche area exponent behaves when the network 
topology depends on sandpile dynamics, i.e. when there 
are mutual interactions between the network structure 
and network dynamics? 

In |llj | Bianconi and Marsili proposed a simple model 
in which the network reorganizes its structure as a con- 
sequence of avalanches of rewiring processes. The only 
parameter of the model which influences the rewiring and 
in consequence the network structure is a type of prob- 
ability that a chosen node becomes unstable and has to 
be rewired. Choosing this probability as a power law 
one can set the system in a critical state and force the 
network to take a power law degree distribution. 

In the present paper, instead of forcing the network to 
stay in a critical state, we allow the system to naturally 
evolve towards the critical region. In our model: 

1. the degree distribution of the considered net- 
works changes due to the distribution of sandpile 
avalanches on this network and 

2. the avalanche size distribution changes because the 
network structure evolves. 

These two mechanisms influence each other and lead to 
the equilibrium point in which the shapes of avalanche 
distribution and degree distribution become identical. In 
the second part of the paper we propose how the pre- 
sented phenomenon can be applied to modelling of evo- 
lution of scientific collaboration. 

In order to complete the rules of our model, apart from 
the rules of sandpile model recapitulated above, we define 
the rewiring process in the following way: each end of a 
link has been assigned a value specifying the time when 
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FIG. 1: Node degree distributions of rewired network in three 
time moments: t sta rt = 0, t mid = 2500 and t end — 10000. 



it was rewired for the last time. After an avalanche of 
area A, the number of A 'oldest' ends of links are rewired 
to the node which has triggered the avalanche off. 

In our studies all networks have: (k) = 4, number of 
nodes N = 10 5 and / = 1(T 4 . We start our simulation 
with Erdos-Renyi random network (which corresponds 
to 7 = oo in static SF networks). Time unit that has 
been used was simply one avalanche and we have car- 
ried out our simulations for t = 10000 steps. Fig. 
presents three snapshots of the node degree distribution 
in three different moments: t s tart — 0, t m j<j = 2500 and 
tend = 10000. Comparing the snapshots, one can see that 
the network reorganizes itself from Poissonian, through a 
mixture of Poissonian and scale- free and finally settles on 
a pure scale-free degree distribution with the well estab- 
lished characteristic exponent 7 = 2. At the same time 
the characteristic exponent of the avalanche area distri- 
bution increases from r = 1.5, which is the known result 
for ER random networks 01, to r = 2. As one can see 
at Fig. |3 in the course of simulation the two exponents 
converge to the common equilibrium value r = 7 = 2. 

Fig. |3| presents the convergence process in a more de- 
tailed way. We define there a new parameter 7(f) that 
in some sense may be understood as the characteristic 
exponent of fat-tailed degree distributions and may be 
compared to 7. j(t) is simply obtained from the sec- 
ond moment of the degree distribution that is know from 
simulations. Given (fc 2 ), we numerically solve the below 
equation for j(t) 
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where 



(7-2)^r(fc- 7 +l,(fc>^) 

Pan(k) = (k) — 



is the known analytic solution of the static model [l7| . 
The new parameter has been introduced because in inter- 
mediate times of the simulation the degree distribution 
does not follow a pure power law. At the mentioned fig- 
ure one can see that the value of the exponent 7 (open 
triangles) decreases from 00 to 2.1. Simultaneously, the 
parameter t characterizing the sandpile dynamics (solid 
squares) increases from 1.5 and finally settles at the equi- 
librium value of 2.1. The values of r fairly good agree 
with the values of Ttheor (solid line) calculated from the 
relation derived for static SF networks by Goh et al. 
[l4| . The last observation let us suspect that during sim- 
ulation the system moves along the trajectory given by 
the formula (see Fig. 0} 
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and respectively the equilibrium point may be calculated 
from the above equation when one assumes 7 = r. 

The exponents 7 = r = 2 characterizing the final crit- 
ical state seem to be robust against different threshold 
assignment strategies in sandpile dynamics. In order to 
support the last statement let us mention two papers 
[H E3 in which we have found probable symptoms of 
such a universality. In |l5j . a class of sandpile mod- 
els is studied. In this class the threshold height of a 
node is set as k l ~ v , where < ri < 1 is a parameter 
of the class. The avalanche size exponent is received 
as t = (7 — 277)/ (7 — 1 — rj). If one assumes that due 
to self-organization and rewiring process in our model 
t = 7, then one obtains r = 7 = 2 independently on 
77. The second example is a model of rapid rearrange- 
ments in the network of the magnetic field flows in the 
Sun corona [lj. The authors show that the avalanches 
of link reconnections and scale free structure of the con- 
sidered network co-generate each other. They also show 
that for the equilibrium the degree distribution exponent 
7 = 2. Unfortunately, they do not present reconnection 
distribution exponent which corresponds to r. 

In fact, the precise value of the fixed point is a bit 
larger than 2, about 2.1. It can result from the finite size 
effects or from the fact, that in the vicinity of 7 = 2 the 
considered networks are highly correlated but the way 
we perform rewiring includes a small random contribu- 
tion (it is known that in this range of parameters the 
network should be correlated disassortatively, so instead 
of rewiring the oldest end of the link we should perhaps 
rewire the least disassorative link). 

In the last part of the paper we would like to show 
how the observed phenomenon can be used to construct 
a simple model of evolution of scientific connections. 

In the model each node represents a scientist. A link 
between two scientists represents the fact that they can 
exchange new ideas and draw inspiration from each other. 
A scientist has its own potential (hidden variable) which 
describes his/her ability to produce an interesting paper. 
The potential is a non decreasing function of time (we 
rather collect ideas, do not lose them). If the potential 
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FIG. 2: Distributions of avalanche area and node degree in 
time t en d- Data are logarithmically binned. Lines are linearly 
fitted with the values indicated at the figure. 
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FIG. 3: Process of equilibration of exponents 7 and r. Solid 
line presents theoretical Ttheor obtained from 7 and eq. 
Inset: second moment of the degree distribution in time. 



crosses a critical level then an interesting publication is 
born and the potential is reduced to zero i.e. the scientist 
has to collect new ideas from the beginning. From time to 
time an avalanche of such interesting publications occurs. 
All those publications use some ideas and draw inspira- 
tion from the first one in the avalanche and of course cite 
it. After the avalanche the scientist who triggered the 
avalanche off attracts attention of other scientists. They 
will try to establish a closer collaboration with such a 
famous person. They will be glad to present to him/her 
their best ideas and they will follow his/her next papers 
because they expect that it is worth to read them. The 
larger avalanche the scientist caused the larger number 
of people will be interested in his/her publications. Be- 
cause it is impossible to observe works of all people every 
scientist reads only the papers of a group of the most im- 
portant persons. That is why after finding a new person 
who become interesting, the scientist has to resign from 
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FIG. 4: Dependence of avalanche area exponent r on parame- 
ter 7 of generated scale- free network. The arrow shows direc- 
tion of a movement in space of parameters during the process 
of equilibration. Measurements are done in equal time steps 
At — 500 and marked as open circles. The black dot depicts 
the fixed point of the process. Dashed line presents theoretical 
Ttheor obtained from 7 and eq. 



observation of other less interesting one. 

In the model that has been presented above one can 
find mechanisms which can be modeled by the phe- 
nomenon described in the first part of the paper. The 
potential may be considered as a column of grains and 
every grain represents a quantum of idea. To make the 
potential equal for all scientists one can normalize it by 
a scientist's degree. By an avalanche we understand an 
occurrence of some scientific sub-domain (like physics 
of complex networks or self-organized criticality) so the 
timescale of such an avalanche will be expressed in years. 

The above model just gives an idea where the pre- 
sented phenomenon can be found. However, since the 
term 'collaboration' does not only mean common papers 
but every manifestation of exchange of ideas (appeared 
as co-authorship, citing and even common discussions), 
the applicability of the model to real data may suffer a 
number of difficulties. 

To conclude, in this paper we have presented by nu- 
merical simulations how the avalanche dynamics of the 
Bak-Tang-Wiesenfeld sandpile model and the network 
structure may influence each other. Such an interplay be- 
tween dynamics and structure leads to self-organization 
in which the shapes of avalanche distribution and degree 
distribution become identical. We suspect that the value 
of both exponents 7 = t = 2 may be universal for a 
large class of SOC phenomena in which the critical be- 
havior occurs not 'on' the network structure but 'in' the 
structure. We also show how the observed phenomenon 
can be used to construct a simple model of evolution of 
scientific collaborations. 
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